# load required libraries
library(tidyverse)
[30m── [1mAttaching packages[22m ────────────────────────────────── tidyverse 1.2.1 ──[39m
[30m[32m✔[30m [34mggplot2[30m 3.1.1 [32m✔[30m [34mpurrr [30m 0.2.5
[32m✔[30m [34mtibble [30m 2.1.1 [32m✔[30m [34mdplyr [30m 0.7.6
[32m✔[30m [34mtidyr [30m 0.8.1 [32m✔[30m [34mstringr[30m 1.3.1
[32m✔[30m [34mreadr [30m 1.1.1 [32m✔[30m [34mforcats[30m 0.3.0[39m
package ‘ggplot2’ was built under R version 3.5.2package ‘tibble’ was built under R version 3.5.2[30m── [1mConflicts[22m ───────────────────────────────────── tidyverse_conflicts() ──
[31m✖[30m [34mdplyr[30m::[32mfilter()[30m masks [34mstats[30m::filter()
[31m✖[30m [34mdplyr[30m::[32mlag()[30m masks [34mstats[30m::lag()[39m
library(langcog) # source: https://github.com/langcog/langcog-package
Attaching package: ‘langcog’
The following object is masked from ‘package:base’:
scale
library(psych)
Attaching package: ‘psych’
The following objects are masked from ‘package:ggplot2’:
%+%, alpha
library(lme4)
Loading required package: Matrix
Attaching package: ‘Matrix’
The following object is masked from ‘package:tidyr’:
expand
library(kableExtra)
# set theme for ggplots
theme_set(theme_bw())
chosen_rot <- "oblimin"
# run source code (extra home-made functions)
source("./scripts/max_factors_efa.R")
source("./scripts/plot_fun.R")
source("./scripts/reten_fun.R")
source("./scripts/data_prep.R")
NAs introduced by coercionattributes are not identical across measure variables;
they will be droppedJoining, by = "question_qualtrics"
Joining, by = "question"
Joining, by = c("ResponseId", "attn_free_coded")
Column `ResponseId` joining character vector and factor, coercing into character vectorColumn `attn_free_coded` joining character vector and factor, coercing into character vectorJoining, by = "ResponseId"
Joining, by = "ResponseId"
NAs introduced by coercionJoining, by = "capacity"
Column `capacity` joining character vector and factor, coercing into character vectorJoining, by = "capacity"
Column `capacity` joining character vector and factor, coercing into character vectorJoining, by = "capacity"
Column `capacity` joining character vector and factor, coercing into character vectorJoining, by = "capacity"
Column `capacity` joining character vector and factor, coercing into character vector
“Baby Mental Life: Study 3” was conducted on MTurk on 2019-04-24.
Our planned sample was 300 participants; based on Study 2, we initially recruited 352 participants. After filtering out participants who failed at least one of our attention checks, we ended up retaining 279 participants (retention rate: 79.3%), but at this point we have not supplemented with additional recruiting. At each stage, we recruited women and men through separate studies, in hopes of acquiring a roughly equal split between genders.
In the end, we ended up with a sample of 279 participants who passed our attention checks, 249 of whom came from unique GPS coordinates (89.2%).
For this first pass, these data include participants where there is another participant with an identical set of GPS coordinates as recorded by Qualtrics.
Each participant assessed children’s mental capacities at 13 target ages between the ages of 0 and 5 years. For each target, they rated 8 mental capacities on a scale from 0 (not at all capable) to 100 (completely capable). In contrast to Study 2, participants completed all 13 ratings of a particular mental capacity on a single “trial” (rather than completing all 8 ratings of a particular target age on a single trial). In addition, participants answered questions about the “developmental factors” that might contribute to development in this domain.
For more details about the study, see our preregistration here.
# load in EFA results from study 1
efa_S1 <- readRDS("../study 1/s1_efa.rds")
heatmap_fun(efa_S1) +
labs(title = paste0("STUDY 1 Parallel Analysis (rotation: ", chosen_rot, ")"),
subtitle = "'% var.' indicates the amount of shared variance explained (total = 100%)")
Joining, by = "capacity"
Joining, by = "factor"
ggplot(d_cap_rating %>%
arrange(ResponseId, domain, capacity, target_ord) %>%
mutate(domain = recode_factor(
domain,
"BOD" = "Bodily sensations",
"NEG" = "Negative emotions",
"POS" = "Social abilities & positive emotions",
"COG" = "Cognition & control",
.default = NA_character_),
capacity = gsub("_", " ", capacity)),
aes(x = target_num, y = response, group = ResponseId,
color = domain)) +
facet_wrap(~ domain ~ capacity, ncol = 4) +
geom_path(alpha = 0.1) +
geom_smooth(aes(group = capacity),
method = "loess",
color = "black") +
scale_x_continuous(breaks = as.numeric(levels(factor(d_cap_rating$target_num))),
labels = levels(d_cap_rating$target_ord)) +
theme(legend.position = "none",
axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
labs(title = "Perceptions of development, by domain",
subtitle = "Black line = loess smoothing",
x = "Target age (numeric)",
y = "Response (0 = not at all capable, 100 = completely capable)")
ggplot(d_cap_rating %>%
arrange(ResponseId, domain, capacity, target_ord) %>%
mutate(domain = recode_factor(
domain,
"BOD" = "Bodily sensations",
"NEG" = "Negative emotions",
"POS" = "Social abilities & positive emotions",
"COG" = "Cognition & control",
.default = NA_character_),
capacity = gsub("_", " ", capacity)),
aes(x = sqrt(target_num), y = response, group = ResponseId,
color = domain)) +
facet_wrap(~ domain ~ capacity, ncol = 4) +
geom_path(alpha = 0.1) +
geom_smooth(aes(group = capacity),
method = "loess",
color = "black") +
scale_x_continuous(breaks = sqrt(as.numeric(levels(factor(d_cap_rating$target_num)))),
labels = levels(d_cap_rating$target_ord)) +
theme(legend.position = "none",
axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
labs(title = "Perceptions of development, by domain",
subtitle = "Black line = loess smoothing",
x = "Target age (numeric, square-root transformed)",
y = "Response (0 = not at all capable, 100 = completely capable)")
ggplot(d_cap_rating %>%
arrange(ResponseId, domain, capacity, target_ord) %>%
mutate(domain = recode_factor(
domain,
"BOD" = "Bodily sensations",
"NEG" = "Negative emotions",
"POS" = "Social abilities & positive emotions",
"COG" = "Cognition & control",
.default = NA_character_),
capacity = gsub("_", " ", capacity)),
aes(x = target_ord, y = response, group = ResponseId,
color = domain)) +
facet_wrap(~ domain ~ capacity, ncol = 4) +
geom_path(alpha = 0.1) +
geom_smooth(aes(group = capacity),
method = "loess",
color = "black") +
theme(legend.position = "none",
axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
labs(title = "Perceptions of development, by domain",
subtitle = "Black line = loess smoothing",
x = "Target age (ordinal)",
y = "Response (0 = not at all capable, 100 = completely capable)")
contrasts_domain_eff_noCOG <- cbind(BOD = c(1, -1, 0, 0),
NEG = c(0, -1, 1, 0),
POS = c(0, -1, 0, 1))
contrasts_domain_eff_noPOS <- cbind(BOD = c(1, 0, 0, -1),
NEG = c(0, 0, 1, -1),
COG = c(0, 1, 0, -1))
contrasts_domain_eff_noNEG <- cbind(BOD = c(1, 0, -1, 0),
POS = c(0, 0, -1, 1),
COG = c(0, 1, -1, 0))
contrasts_domain_eff_noBOD <- cbind(POS = c(-1, 0, 0, 1),
NEG = c(-1, 0, 1, 0),
COG = c(-1, 1, 0, 0))
# r1_noBOD <- lmer(response ~ target_num * domain
# + (target_num + domain | ResponseId)
# + (target_num | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12),
# contrasts = list(domain = contrasts_domain_eff_noBOD))
# saveRDS(r1_noBOD, "./models/r1_noBOD.RDS")
r1_noBOD <- readRDS("./models/r1_noBOD.RDS")
summary(r1_noBOD, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula:
response ~ target_num * domain + (target_num + domain | ResponseId) +
(target_num | capacity)
Data: d_cap_rating %>% mutate(target_num = target_num/12)
REML criterion at convergence: 228041.6
Scaled residuals:
Min 1Q Median 3Q Max
-5.2123 -0.4827 0.0427 0.5204 4.6348
Random effects:
Groups Name Variance Std.Dev. Corr
ResponseId (Intercept) 195.361 13.977
target_num 9.655 3.107 -0.69
domainPOS 91.876 9.585 0.34 -0.29
domainNEG 196.219 14.008 0.57 -0.41 -0.24
domainCOG 175.285 13.240 -0.22 0.30 -0.24 -0.58
capacity (Intercept) 193.647 13.916
target_num 6.995 2.645 -0.86
Residual 345.942 18.600
Number of obs: 25896, groups: ResponseId, 249; capacity, 8
Fixed effects:
Estimate Std. Error t value
(Intercept) 55.2190 5.0015 11.041
target_num 8.1121 0.9584 8.465
domainPOS -3.8395 8.5475 -0.449
domainNEG 9.0385 8.5719 1.054
domainCOG -44.3093 8.5670 -5.172
target_num:domainPOS 3.8509 1.6245 2.370
target_num:domainNEG -2.3475 1.6245 -1.445
target_num:domainCOG 5.4235 1.6245 3.339
# r2_noBOD <- lmer(response ~ poly(target_num, 3) * domain
# + (target_num + domain | ResponseId)
# + (target_num | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12),
# contrasts = list(domain = contrasts_domain_eff_noBOD))
# saveRDS(r2_noBOD, "./models/r2_noBOD.RDS")
r2_noBOD <- readRDS("./models/r2_noBOD.RDS")
summary(r2_noBOD, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula: response ~ poly(target_num, 3) * domain + (target_num + domain |
ResponseId) + (target_num | capacity)
Data: d_cap_rating %>% mutate(target_num = target_num/12)
REML criterion at convergence: 223981
Scaled residuals:
Min 1Q Median 3Q Max
-5.2929 -0.4858 0.0105 0.5052 5.3522
Random effects:
Groups Name Variance Std.Dev. Corr
ResponseId (Intercept) 89.405 9.455
target_num 9.853 3.139 -0.56
domainCOG 374.649 19.356 -0.47 -0.02
domainNEG 435.336 20.865 0.23 -0.48 0.27
domainPOS 277.020 16.644 0.04 -0.43 0.50 0.51
capacity (Intercept) 193.770 13.920
target_num 7.004 2.646 -0.86
Residual 294.790 17.169
Number of obs: 25896, groups: ResponseId, 249; capacity, 8
Fixed effects:
Estimate Std. Error t value
(Intercept) 66.666 3.914 17.031
poly(target_num, 3)1 2058.398 243.332 8.459
poly(target_num, 3)2 -844.244 17.169 -49.171
poly(target_num, 3)3 407.030 17.169 23.707
domainPOS 1.595 6.691 0.238
domainNEG 5.726 6.723 0.852
domainCOG -36.656 6.716 -5.458
poly(target_num, 3)1:domainPOS 977.144 412.297 2.370
poly(target_num, 3)2:domainPOS -862.372 29.738 -28.999
poly(target_num, 3)3:domainPOS 507.310 29.738 17.059
poly(target_num, 3)1:domainNEG -595.664 412.297 -1.445
poly(target_num, 3)2:domainNEG 173.963 29.738 5.850
poly(target_num, 3)3:domainNEG -92.256 29.738 -3.102
poly(target_num, 3)1:domainCOG 1376.179 412.297 3.338
poly(target_num, 3)2:domainCOG 53.207 29.738 1.789
poly(target_num, 3)3:domainCOG -157.616 29.738 -5.300
Correlation matrix not shown by default, as p = 16 > 12.
Use print(x, correlation=TRUE) or
vcov(x) if you need it
# r1_noNEG <- lmer(response ~ target_num * domain
# + (target_num + domain | ResponseId)
# + (target_num | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12),
# contrasts = list(domain = contrasts_domain_eff_noNEG))
# saveRDS(r1_noNEG, "./models/r1_noNEG.RDS")
r1_noNEG <- readRDS("./models/r1_noNEG.RDS")
summary(r1_noNEG, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula:
response ~ target_num * domain + (target_num + domain | ResponseId) +
(target_num | capacity)
Data: d_cap_rating %>% mutate(target_num = target_num/12)
REML criterion at convergence: 228041.6
Scaled residuals:
Min 1Q Median 3Q Max
-5.2123 -0.4827 0.0427 0.5204 4.6348
Random effects:
Groups Name Variance Std.Dev. Corr
ResponseId (Intercept) 195.361 13.977
target_num 9.655 3.107 -0.69
domainPOS 91.876 9.585 0.34 -0.29
domainNEG 196.220 14.008 0.57 -0.41 -0.24
domainCOG 175.286 13.240 -0.22 0.30 -0.24 -0.58
capacity (Intercept) 193.749 13.919
target_num 6.997 2.645 -0.86
Residual 345.942 18.600
Number of obs: 25896, groups: ResponseId, 249; capacity, 8
Fixed effects:
Estimate Std. Error t value
(Intercept) 55.2190 5.0027 11.038
target_num 8.1121 0.9585 8.463
domainBOD 39.1102 8.5568 4.571
domainPOS -3.8395 8.5497 -0.449
domainCOG -44.3093 8.5692 -5.171
target_num:domainBOD -6.9269 1.6248 -4.263
target_num:domainPOS 3.8509 1.6248 2.370
target_num:domainCOG 5.4235 1.6248 3.338
# r2_noNEG <- lmer(response ~ poly(target_num, 3) * domain
# + (target_num + domain | ResponseId)
# + (target_num | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12),
# contrasts = list(domain = contrasts_domain_eff_noNEG))
# saveRDS(r2_noNEG, "./models/r2_noNEG.RDS")
r2_noNEG <- readRDS("./models/r2_noNEG.RDS")
summary(r2_noNEG, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula: response ~ poly(target_num, 3) * domain + (target_num + domain |
ResponseId) + (target_num | capacity)
Data: d_cap_rating %>% mutate(target_num = target_num/12)
REML criterion at convergence: 223981
Scaled residuals:
Min 1Q Median 3Q Max
-5.2929 -0.4858 0.0105 0.5052 5.3522
Random effects:
Groups Name Variance Std.Dev. Corr
ResponseId (Intercept) 89.406 9.455
target_num 9.853 3.139 -0.56
domainCOG 374.656 19.356 -0.47 -0.02
domainNEG 435.338 20.865 0.23 -0.48 0.27
domainPOS 277.024 16.644 0.04 -0.43 0.50 0.51
capacity (Intercept) 193.741 13.919
target_num 7.003 2.646 -0.86
Residual 294.790 17.169
Number of obs: 25896, groups: ResponseId, 249; capacity, 8
Fixed effects:
Estimate Std. Error t value
(Intercept) 66.666 3.914 17.032
poly(target_num, 3)1 2058.398 243.316 8.460
poly(target_num, 3)2 -844.244 17.169 -49.171
poly(target_num, 3)3 407.030 17.169 23.707
domainBOD 29.336 6.700 4.378
domainPOS 1.595 6.691 0.238
domainCOG -36.656 6.716 -5.458
poly(target_num, 3)1:domainBOD -1757.660 412.269 -4.263
poly(target_num, 3)2:domainBOD 635.202 29.738 21.360
poly(target_num, 3)3:domainBOD -257.438 29.738 -8.657
poly(target_num, 3)1:domainPOS 977.144 412.269 2.370
poly(target_num, 3)2:domainPOS -862.372 29.738 -28.999
poly(target_num, 3)3:domainPOS 507.310 29.738 17.059
poly(target_num, 3)1:domainCOG 1376.179 412.269 3.338
poly(target_num, 3)2:domainCOG 53.207 29.738 1.789
poly(target_num, 3)3:domainCOG -157.616 29.738 -5.300
Correlation matrix not shown by default, as p = 16 > 12.
Use print(x, correlation=TRUE) or
vcov(x) if you need it
# r1_noCOG <- lmer(response ~ target_num * domain
# + (target_num + domain | ResponseId)
# + (target_num | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12),
# contrasts = list(domain = contrasts_domain_eff_noCOG))
# saveRDS(r1_noCOG, "./models/r1_noCOG.RDS")
r1_noCOG <- readRDS("./models/r1_noCOG.RDS")
summary(r1_noCOG, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula:
response ~ target_num * domain + (target_num + domain | ResponseId) +
(target_num | capacity)
Data: d_cap_rating %>% mutate(target_num = target_num/12)
REML criterion at convergence: 228041.6
Scaled residuals:
Min 1Q Median 3Q Max
-5.2123 -0.4827 0.0427 0.5204 4.6348
Random effects:
Groups Name Variance Std.Dev. Corr
ResponseId (Intercept) 195.360 13.977
target_num 9.655 3.107 -0.69
domainPOS 91.877 9.585 0.34 -0.29
domainNEG 196.219 14.008 0.57 -0.41 -0.24
domainCOG 175.284 13.239 -0.22 0.30 -0.24 -0.58
capacity (Intercept) 193.602 13.914
target_num 6.993 2.644 -0.86
Residual 345.942 18.600
Number of obs: 25896, groups: ResponseId, 249; capacity, 8
Fixed effects:
Estimate Std. Error t value
(Intercept) 55.2190 5.0009 11.042
target_num 8.1121 0.9583 8.465
domainBOD 39.1102 8.5536 4.572
domainNEG 9.0385 8.5709 1.055
domainPOS -3.8395 8.5464 -0.449
target_num:domainBOD -6.9269 1.6244 -4.264
target_num:domainNEG -2.3475 1.6244 -1.445
target_num:domainPOS 3.8509 1.6244 2.371
# r2_noCOG <- lmer(response ~ poly(target_num, 3) * domain
# + (target_num + domain | ResponseId)
# + (target_num | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12),
# contrasts = list(domain = contrasts_domain_eff_noCOG))
# saveRDS(r2_noCOG, "./models/r2_noCOG.RDS")
r2_noCOG <- readRDS("./models/r2_noCOG.RDS")
summary(r2_noCOG, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula: response ~ poly(target_num, 3) * domain + (target_num + domain |
ResponseId) + (target_num | capacity)
Data: d_cap_rating %>% mutate(target_num = target_num/12)
REML criterion at convergence: 223981
Scaled residuals:
Min 1Q Median 3Q Max
-5.2929 -0.4858 0.0105 0.5052 5.3522
Random effects:
Groups Name Variance Std.Dev. Corr
ResponseId (Intercept) 89.405 9.455
target_num 9.853 3.139 -0.56
domainCOG 374.650 19.356 -0.47 -0.02
domainNEG 435.333 20.865 0.23 -0.48 0.27
domainPOS 277.018 16.644 0.04 -0.43 0.50 0.51
capacity (Intercept) 193.787 13.921
target_num 7.004 2.646 -0.86
Residual 294.790 17.169
Number of obs: 25896, groups: ResponseId, 249; capacity, 8
Fixed effects:
Estimate Std. Error t value
(Intercept) 66.666 3.915 17.030
poly(target_num, 3)1 2058.398 243.334 8.459
poly(target_num, 3)2 -844.244 17.169 -49.171
poly(target_num, 3)3 407.030 17.169 23.707
domainBOD 29.336 6.701 4.378
domainNEG 5.726 6.723 0.852
domainPOS 1.595 6.692 0.238
poly(target_num, 3)1:domainBOD -1757.660 412.300 -4.263
poly(target_num, 3)2:domainBOD 635.202 29.738 21.360
poly(target_num, 3)3:domainBOD -257.438 29.738 -8.657
poly(target_num, 3)1:domainNEG -595.664 412.300 -1.445
poly(target_num, 3)2:domainNEG 173.963 29.738 5.850
poly(target_num, 3)3:domainNEG -92.256 29.738 -3.102
poly(target_num, 3)1:domainPOS 977.144 412.300 2.370
poly(target_num, 3)2:domainPOS -862.372 29.738 -28.999
poly(target_num, 3)3:domainPOS 507.310 29.738 17.059
Correlation matrix not shown by default, as p = 16 > 12.
Use print(x, correlation=TRUE) or
vcov(x) if you need it
NOTE: All of these models fail to meet our criterion of random effects being correlation < 0.90. Because of this, I have not included in this print-out and I have not yet modeled non-linear effects after square-root transformation. The next step will be to try simpler random effects models (with fewer random slopes).
NOTE: All of the preregistered models fail to meet our criterion of random effects being correlation < 0.90. Because of this, I instead used simpler random effects models (with fewer random slopes), as indicated in the preregistration.
# r5a_BOD <- lmer(response ~ target_num
# + (1 | ResponseId)
# + (1 | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12) %>%
# filter(domain == "BOD"))
# saveRDS(r5a_BOD, "./models/r5a_BOD.RDS")
r5a_BOD <- readRDS("./models/r5a_BOD.RDS")
summary(r5a_BOD, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula: response ~ target_num + (1 | ResponseId) + (1 | capacity)
Data:
d_cap_rating %>% mutate(target_num = target_num/12) %>% filter(domain ==
"BOD")
REML criterion at convergence: 55342
Scaled residuals:
Min 1Q Median 3Q Max
-8.0547 -0.1326 0.1032 0.2568 3.8043
Random effects:
Groups Name Variance Std.Dev.
ResponseId (Intercept) 101.283 10.064
capacity (Intercept) 1.223 1.106
Residual 106.243 10.307
Number of obs: 7254, groups: ResponseId, 279; capacity, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) 93.31343 1.00055 93.26
target_num 1.32602 0.07675 17.28
# r6a_BOD <- lmer(response ~ poly(target_num, 3)
# + (1 | ResponseId)
# + (1 | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12) %>%
# filter(domain == "BOD"))
# saveRDS(r6a_BOD, "./models/r6a_BOD.RDS")
r6a_BOD <- readRDS("./models/r6a_BOD.RDS")
summary(r6a_BOD, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula:
response ~ poly(target_num, 3) + (1 | ResponseId) + (1 | capacity)
Data:
d_cap_rating %>% mutate(target_num = target_num/12) %>% filter(domain ==
"BOD")
REML criterion at convergence: 48865.4
Scaled residuals:
Min 1Q Median 3Q Max
-8.1452 -0.2055 0.0208 0.3275 3.8911
Random effects:
Groups Name Variance Std.Dev.
ResponseId (Intercept) 69.69 8.348
capacity (Intercept) 1.08 1.039
Residual 99.44 9.972
Number of obs: 6474, groups: ResponseId, 249; capacity, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) 96.002 0.914 105.038
poly(target_num, 3)1 150.369 9.972 15.079
poly(target_num, 3)2 -104.521 9.972 -10.482
poly(target_num, 3)3 74.796 9.972 7.501
# r5a_NEG <- lmer(response ~ target_num
# + (1 | ResponseId)
# + (1 | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12) %>%
# filter(domain == "NEG"))
# saveRDS(r5a_NEG, "./models/r5a_NEG.RDS")
r5a_NEG <- readRDS("./models/r5a_NEG.RDS")
summary(r5a_NEG, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula: response ~ target_num + (1 | ResponseId) + (1 | capacity)
Data:
d_cap_rating %>% mutate(target_num = target_num/12) %>% filter(domain ==
"NEG")
REML criterion at convergence: 67536.4
Scaled residuals:
Min 1Q Median 3Q Max
-3.8357 -0.5649 0.0850 0.7065 3.1869
Random effects:
Groups Name Variance Std.Dev.
ResponseId (Intercept) 485.5 22.03
capacity (Intercept) 233.9 15.30
Residual 573.1 23.94
Number of obs: 7254, groups: ResponseId, 279; capacity, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) 63.4263 10.9020 5.818
target_num 5.9054 0.1783 33.128
# r6a_NEG <- lmer(response ~ poly(target_num, 3)
# + (1 | ResponseId)
# + (1 | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12) %>%
# filter(domain == "NEG"))
# saveRDS(r6a_NEG, "./models/r6a_NEG.RDS")
r6a_NEG <- readRDS("./models/r6a_NEG.RDS")
summary(r6a_NEG, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula:
response ~ poly(target_num, 3) + (1 | ResponseId) + (1 | capacity)
Data:
d_cap_rating %>% mutate(target_num = target_num/12) %>% filter(domain ==
"NEG")
REML criterion at convergence: 59981.5
Scaled residuals:
Min 1Q Median 3Q Max
-3.6632 -0.6339 0.0748 0.6934 3.4866
Random effects:
Groups Name Variance Std.Dev.
ResponseId (Intercept) 491.7 22.18
capacity (Intercept) 232.9 15.26
Residual 549.0 23.43
Number of obs: 6474, groups: ResponseId, 249; capacity, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) 72.39 10.89 6.650
poly(target_num, 3)1 731.37 23.43 31.214
poly(target_num, 3)2 -335.14 23.43 -14.304
poly(target_num, 3)3 157.39 23.43 6.717
# r5a_COG <- lmer(response ~ target_num
# + (1 | ResponseId)
# + (1 | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12) %>%
# filter(domain == "COG"))
# saveRDS(r5a_COG, "./models/r5a_COG.RDS")
r5a_COG <- readRDS("./models/r5a_COG.RDS")
summary(r5a_COG, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula: response ~ target_num + (1 | ResponseId) + (1 | capacity)
Data:
d_cap_rating %>% mutate(target_num = target_num/12) %>% filter(domain ==
"COG")
REML criterion at convergence: 62229.9
Scaled residuals:
Min 1Q Median 3Q Max
-5.1523 -0.6120 -0.0335 0.5330 4.1114
Random effects:
Groups Name Variance Std.Dev.
ResponseId (Intercept) 248.02 15.749
capacity (Intercept) 25.34 5.033
Residual 275.14 16.587
Number of obs: 7254, groups: ResponseId, 279; capacity, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) 10.7635 3.6912 2.916
target_num 13.4807 0.1235 109.145
# r6a_COG <- lmer(response ~ poly(target_num, 3)
# + (1 | ResponseId)
# + (1 | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12) %>%
# filter(domain == "COG"))
# saveRDS(r6a_COG, "./models/r6a_COG.RDS")
r6a_COG <- readRDS("./models/r6a_COG.RDS")
summary(r6a_COG, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula:
response ~ poly(target_num, 3) + (1 | ResponseId) + (1 | capacity)
Data:
d_cap_rating %>% mutate(target_num = target_num/12) %>% filter(domain ==
"COG")
REML criterion at convergence: 54889.4
Scaled residuals:
Min 1Q Median 3Q Max
-4.9686 -0.6187 -0.0333 0.5887 4.1843
Random effects:
Groups Name Variance Std.Dev.
ResponseId (Intercept) 261.07 16.158
capacity (Intercept) 29.69 5.449
Residual 248.48 15.763
Number of obs: 6474, groups: ResponseId, 249; capacity, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) 30.010 3.991 7.519
poly(target_num, 3)1 1717.288 15.763 108.942
poly(target_num, 3)2 -395.519 15.763 -25.091
poly(target_num, 3)3 124.707 15.763 7.911
# r7a_BOD <- lmer(response ~ sqrt(target_num)
# + (1 | ResponseId)
# + (1 | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12) %>%
# filter(domain == "BOD"))
# saveRDS(r7a_BOD, "./models/r7a_BOD.RDS")
r7a_BOD <- readRDS("./models/r7a_BOD.RDS")
summary(r7a_BOD, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula: response ~ sqrt(target_num) + (1 | ResponseId) + (1 | capacity)
Data:
d_cap_rating %>% mutate(target_num = target_num/12) %>% filter(domain ==
"BOD")
REML criterion at convergence: 55196.4
Scaled residuals:
Min 1Q Median 3Q Max
-7.9734 -0.2182 0.0823 0.2948 4.0115
Random effects:
Groups Name Variance Std.Dev.
ResponseId (Intercept) 101.366 10.068
capacity (Intercept) 1.224 1.106
Residual 104.071 10.202
Number of obs: 7254, groups: ResponseId, 279; capacity, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) 91.6317 1.0088 90.83
sqrt(target_num) 3.6774 0.1733 21.22
# r8a_BOD <- lmer(response ~ poly(sqrt(target_num), 3)
# + (1 | ResponseId)
# + (1 | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12) %>%
# filter(domain == "BOD"))
# saveRDS(r8a_BOD, "./models/r8a_BOD.RDS")
r8a_BOD <- readRDS("./models/r8a_BOD.RDS")
summary(r8a_BOD, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula: response ~ poly(sqrt(target_num), 3) + (1 | ResponseId) + (1 |
capacity)
Data:
d_cap_rating %>% mutate(target_num = target_num/12) %>% filter(domain ==
"BOD")
REML criterion at convergence: 48802.4
Scaled residuals:
Min 1Q Median 3Q Max
-8.0087 -0.1792 -0.0219 0.3389 4.1311
Random effects:
Groups Name Variance Std.Dev.
ResponseId (Intercept) 69.73 8.351
capacity (Intercept) 1.08 1.039
Residual 98.43 9.921
Number of obs: 6474, groups: ResponseId, 249; capacity, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) 96.0017 0.9141 105.028
poly(sqrt(target_num), 3)1 185.1919 9.9214 18.666
poly(sqrt(target_num), 3)2 -98.7964 9.9214 -9.958
poly(sqrt(target_num), 3)3 36.1618 9.9214 3.645
# r7a_NEG <- lmer(response ~ sqrt(target_num)
# + (1 | ResponseId)
# + (1 | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12) %>%
# filter(domain == "NEG"))
# saveRDS(r7a_NEG, "./models/r7a_NEG.RDS")
r7a_NEG <- readRDS("./models/r7a_NEG.RDS")
summary(r7a_NEG, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula: response ~ sqrt(target_num) + (1 | ResponseId) + (1 | capacity)
Data:
d_cap_rating %>% mutate(target_num = target_num/12) %>% filter(domain ==
"NEG")
REML criterion at convergence: 67281.1
Scaled residuals:
Min 1Q Median 3Q Max
-3.8115 -0.5996 0.0806 0.6890 3.5062
Random effects:
Groups Name Variance Std.Dev.
ResponseId (Intercept) 486.3 22.05
capacity (Intercept) 233.9 15.30
Residual 552.6 23.51
Number of obs: 7254, groups: ResponseId, 279; capacity, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) 57.3409 10.9058 5.258
sqrt(target_num) 14.9239 0.3994 37.370
# r8a_NEG <- lmer(response ~ poly(sqrt(target_num), 3)
# + (1 | ResponseId)
# + (1 | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12) %>%
# filter(domain == "NEG"))
# saveRDS(r8a_NEG, "./models/r8a_NEG.RDS")
r8a_NEG <- readRDS("./models/r8a_NEG.RDS")
summary(r8a_NEG, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula: response ~ poly(sqrt(target_num), 3) + (1 | ResponseId) + (1 |
capacity)
Data:
d_cap_rating %>% mutate(target_num = target_num/12) %>% filter(domain ==
"NEG")
REML criterion at convergence: 59962.5
Scaled residuals:
Min 1Q Median 3Q Max
-3.6297 -0.6396 0.0635 0.6857 3.6271
Random effects:
Groups Name Variance Std.Dev.
ResponseId (Intercept) 491.8 22.18
capacity (Intercept) 232.9 15.26
Residual 547.3 23.39
Number of obs: 6474, groups: ResponseId, 249; capacity, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) 72.39 10.89 6.650
poly(sqrt(target_num), 3)1 809.13 23.39 34.586
poly(sqrt(target_num), 3)2 -164.47 23.39 -7.030
poly(sqrt(target_num), 3)3 -25.42 23.39 -1.087
# r7a_COG <- lmer(response ~ sqrt(target_num)
# + (1 | ResponseId)
# + (1 | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12) %>%
# filter(domain == "COG"))
# saveRDS(r7a_COG, "./models/r7a_COG.RDS")
r7a_COG <- readRDS("./models/r7a_COG.RDS")
summary(r7a_COG, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula: response ~ sqrt(target_num) + (1 | ResponseId) + (1 | capacity)
Data:
d_cap_rating %>% mutate(target_num = target_num/12) %>% filter(domain ==
"COG")
REML criterion at convergence: 61697.6
Scaled residuals:
Min 1Q Median 3Q Max
-4.8211 -0.6330 -0.0502 0.5957 3.9643
Random effects:
Groups Name Variance Std.Dev.
ResponseId (Intercept) 248.79 15.773
capacity (Intercept) 25.34 5.034
Residual 254.98 15.968
Number of obs: 7254, groups: ResponseId, 279; capacity, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) -0.5581 3.6968 -0.151
sqrt(target_num) 31.4079 0.2713 115.783
# r8a_COG <- lmer(response ~ poly(sqrt(target_num), 3)
# + (1 | ResponseId)
# + (1 | capacity),
# d_cap_rating %>%
# mutate(target_num = target_num/12) %>%
# filter(domain == "COG"))
# saveRDS(r8a_COG, "./models/r8a_COG.RDS")
r8a_COG <- readRDS("./models/r8a_COG.RDS")
summary(r8a_COG, corr = F)
Linear mixed model fit by REML ['lmerMod']
Formula: response ~ poly(sqrt(target_num), 3) + (1 | ResponseId) + (1 |
capacity)
Data:
d_cap_rating %>% mutate(target_num = target_num/12) %>% filter(domain ==
"COG")
REML criterion at convergence: 54905.8
Scaled residuals:
Min 1Q Median 3Q Max
-4.8884 -0.6120 -0.0278 0.5740 4.1611
Random effects:
Groups Name Variance Std.Dev.
ResponseId (Intercept) 261.04 16.157
capacity (Intercept) 29.69 5.449
Residual 249.14 15.784
Number of obs: 6474, groups: ResponseId, 249; capacity, 2
Fixed effects:
Estimate Std. Error t value
(Intercept) 30.010 3.991 7.519
poly(sqrt(target_num), 3)1 1756.660 15.784 111.293
poly(sqrt(target_num), 3)2 108.002 15.784 6.842
poly(sqrt(target_num), 3)3 -139.527 15.784 -8.840
ggplot(d_dev_factor_rating %>%
mutate(domain = recode_factor(
domain,
"BOD" = "Bodily sensations",
"NEG" = "Negative emotions",
"POS" = "Social abilities & positive emotions",
"COG" = "Cognition & control",
.default = NA_character_)) %>%
mutate_at(vars(capacity, dev_factor),
funs(gsub("_", " ", .))) %>%
mutate(dev_factor = factor(
dev_factor,
levels = gsub("_", " ", levels(d_dev_factor_rating$dev_factor)))),
aes(x = dev_factor, y = response, color = domain)) +
facet_wrap(~ domain ~ capacity, ncol = 4) +
geom_jitter(alpha = 0.025, height = 0.2, width = 0.2) +
geom_pointrange(data = . %>% group_by(domain, capacity, dev_factor) %>%
multi_boot_standard(col = "response", na.rm = T),
aes(y = mean, ymin = ci_lower, ymax = ci_upper),
color = "black", fatten = 4) +
theme(legend.position = "none",
axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
labs(title = "Ratings of developmental factors (Version 1)",
subtitle = "Error bars are bootstrapped 95% CIs",
x = "Developmental factor",
y = "Response (0 = plays no role, 6 = plays a very important role)")
ggplot(d_dev_factor_rating %>%
mutate(domain = recode_factor(
domain,
"BOD" = "Bodily sensations",
"NEG" = "Negative emotions",
"POS" = "Social abilities & positive emotions",
"COG" = "Cognition & control",
.default = NA_character_)) %>%
mutate_at(vars(capacity, dev_factor),
funs(gsub("_", " ", .))) %>%
mutate(dev_factor = factor(
dev_factor,
levels = gsub("_", " ", levels(d_dev_factor_rating$dev_factor)))),
aes(x = reorder(capacity, as.numeric(domain)),
y = response, color = domain)) +
facet_wrap(~ dev_factor, ncol = 5) +
geom_jitter(alpha = 0.025, height = 0.2, width = 0.2) +
geom_pointrange(data = . %>% group_by(domain, capacity, dev_factor) %>%
multi_boot_standard(col = "response", na.rm = T),
aes(y = mean, ymin = ci_lower, ymax = ci_upper),
color = "black", fatten = 4) +
theme(legend.position = "top",
axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
guides(color = guide_legend(override.aes = list(alpha = 1, size = 3))) +
labs(title = "Ratings of developmental factors (Version 2)",
subtitle = "Error bars are bootstrapped 95% CIs",
color = "Domain",
x = "Capacity (by domain)",
y = "Response (0 = plays no role, 6 = plays a very important role)")
ggplot(d_dev_factor_rating %>%
mutate(dev_factor = factor(
gsub("_", " ", dev_factor),
levels = gsub("_", " ", levels(d_dev_factor_rating$dev_factor)))),
aes(x = dev_factor, y = response)) +
geom_jitter(alpha = 0.01, height = 0.3, width = 0.3, color = "blue") +
geom_pointrange(data = . %>% group_by(dev_factor) %>%
multi_boot_standard(col = "response", na.rm = T),
aes(y = mean, ymin = ci_lower, ymax = ci_upper),
color = "black", fatten = 2) +
theme(legend.position = "none",
axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
labs(title = "Ratings of developmental factors (Version 3)",
subtitle = "Error bars are bootstrapped 95% CIs",
x = "Developmental factor",
y = "Response (0 = plays no role, 6 = plays a very important role)")
It could be nice to reduce the number of “devleopmental factors” from 10 down to something more manageable for analysis. Here I explore EFA and hierarchical clustering as dimension reduction methods, and re-plot with these results in mind.
temp <- d_dev_factor_rating %>%
unite(key, ResponseId, domain, capacity) %>%
spread(dev_factor, response) %>%
column_to_rownames("key")
reten_fun(temp, "varimax")
[1] 2
# fa.parallel(temp)
# VSS(temp)
temp_efa <- fa(temp,
# nfactors = 4,
nfactors = reten_fun(temp, "varimax"),
rotate = "varimax")
heatmap_fun(temp_efa) +
guides(fill = guide_colorbar(barheight = 8, barwidth = 0.5))
Joining, by = "capacity"
Joining, by = "factor"
library(ggdendro)
library(dendextend)
---------------------
Welcome to dendextend version 1.8.0
Type citation('dendextend') for how to cite the package.
Type browseVignettes(package = 'dendextend') for the package vignette.
The github page is: https://github.com/talgalili/dendextend/
Suggestions and bug-reports can be submitted at: https://github.com/talgalili/dendextend/issues
Or contact: <tal.galili@gmail.com>
To suppress this message use: suppressPackageStartupMessages(library(dendextend))
---------------------
Attaching package: ‘dendextend’
The following object is masked from ‘package:ggdendro’:
theme_dendro
The following object is masked from ‘package:stats’:
cutree
temp_hclust <- temp %>%
t() %>%
dist() %>%
hclust()
temp_hclust_order <- data.frame(order = as.numeric(temp_hclust$order),
dev_factor = as.character(temp_hclust$labels))
temp_hclust %>%
ggdendrogram(rotate = F) +
theme_minimal() +
theme(#axis.title = element_blank(),
# axis.text.y = element_blank(),
axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1),
panel.grid.major.x = element_blank(),
panel.grid.minor.x = element_blank()) +
labs(title = "Hierarchical agglomerative clustering",
subtitle = "Complete linkage (default for stats::hclust() function)",
x = "Developmental factor", y = "Height")
# as.dendrogram() %>%
# set("labels_col", k = 2) %>%
# plot(horiz = T, xlim = c(-100, 200))
temp2 <- d_dev_factor_rating %>%
mutate(dev_factor_cluster = case_when(
dev_factor %in% c("experiments", "people_teach", "brain_changes",
"observes_objects", "observes_people",
"interacts_people") ~ "external",
dev_factor %in% c("preprogrammed", "womb_experiences", "body_grows",
"senses_improve") ~ "internal",
TRUE ~ NA_character_)) %>%
mutate(dev_factor_cluster = factor(dev_factor_cluster),
response_cent = response - 3) %>%
left_join(temp_hclust_order)
Joining, by = "dev_factor"
Column `dev_factor` joining factors with different levels, coercing to character vector
contrasts_cluster_eff <- cbind(EXT = c(1, -1))
ggplot(temp2 %>%
mutate(domain = recode_factor(
domain,
"BOD" = "Bodily sensations",
"NEG" = "Negative emotions",
"POS" = "Social abilities & positive emotions",
"COG" = "Cognition & control",
.default = NA_character_)) %>%
mutate_at(vars(capacity, dev_factor),
funs(gsub("_", " ", .))),
aes(x = reorder(dev_factor, as.numeric(dev_factor_cluster)),
y = response, color = domain, shape = dev_factor_cluster)) +
facet_wrap(~ domain ~ capacity, ncol = 4) +
geom_jitter(alpha = 0.025, height = 0.2, width = 0.2, show.legend = F) +
geom_pointrange(data = . %>% group_by(domain, capacity,
dev_factor_cluster, dev_factor) %>%
multi_boot_standard(col = "response", na.rm = T),
aes(y = mean, ymin = ci_lower, ymax = ci_upper),
color = "black", fatten = 4) +
theme(legend.position = "top",
axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
labs(title = "Ratings of developmental factors (Version 1)",
subtitle = "Error bars are bootstrapped 95% CIs",
shape = "Type of developmental factor (per EFA/clustering)",
x = "Developmental factor",
y = "Response (0 = plays no role, 6 = plays a very important role)")
ggplot(temp2 %>%
mutate(domain = recode_factor(
domain,
"BOD" = "Bodily sensations",
"NEG" = "Negative emotions",
"POS" = "Social abilities & positive emotions",
"COG" = "Cognition & control",
.default = NA_character_)) %>%
mutate_at(vars(capacity),
funs(gsub("_", " ", .))),
aes(x = reorder(capacity, as.numeric(domain)),
y = response, color = domain)) +
facet_wrap(dev_factor_cluster ~ dev_factor, ncol = 6) +
geom_jitter(alpha = 0.025, height = 0.2, width = 0.2) +
geom_pointrange(data = . %>% group_by(domain, capacity,
dev_factor_cluster, dev_factor) %>%
multi_boot_standard(col = "response", na.rm = T),
aes(y = mean, ymin = ci_lower, ymax = ci_upper),
color = "black", fatten = 4) +
theme(legend.position = "top",
axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
guides(color = guide_legend(override.aes = list(alpha = 1, size = 3))) +
labs(title = "Ratings of developmental factors (Version 3)",
subtitle = "Error bars are bootstrapped 95% CIs",
color = "Domain",
x = "Capacity (by domain)",
y = "Response (0 = plays no role, 6 = plays a very important role)")
r_temp <- lmer(response ~ domain * dev_factor_cluster +
# (domain + dev_factor_cluster | ResponseId) +
# including both random intercepts yields corr > 0.9
(1 + dev_factor_cluster | ResponseId) +
(1 | capacity) + (1 | dev_factor),
temp2,
contrasts = list(domain = contrasts_domain_eff_noBOD,
dev_factor_cluster = contrasts_cluster_eff))
summary(r_temp)
Linear mixed model fit by REML ['lmerMod']
Formula:
response ~ domain * dev_factor_cluster + (1 + dev_factor_cluster |
ResponseId) + (1 | capacity) + (1 | dev_factor)
Data: temp2
REML criterion at convergence: 87685.9
Scaled residuals:
Min 1Q Median 3Q Max
-3.4691 -0.7172 0.0083 0.7176 3.4754
Random effects:
Groups Name Variance Std.Dev. Corr
ResponseId (Intercept) 0.72087 0.8490
dev_factor_clusterinternal 0.31234 0.5589 -0.16
dev_factor (Intercept) 0.46942 0.6851
capacity (Intercept) 0.08958 0.2993
Residual 2.81132 1.6767
Number of obs: 22320, groups:
ResponseId, 279; dev_factor, 10; capacity, 8
Fixed effects:
Estimate Std. Error t value
(Intercept) 3.13801 0.25063 12.520
domainPOS 0.56263 0.18435 3.052
domainNEG -0.22642 0.18435 -1.228
domainCOG 0.22527 0.18435 1.222
dev_factor_clusterEXT 0.35587 0.22206 1.603
domainPOS:dev_factor_clusterEXT 0.34522 0.01984 17.401
domainNEG:dev_factor_clusterEXT -0.10109 0.01984 -5.095
domainCOG:dev_factor_clusterEXT 0.83940 0.01984 42.309
Correlation of Fixed Effects:
(Intr) dmnPOS dmnNEG dmnCOG d__EXT dPOS:_ dNEG:_
domainPOS 0.000
domainNEG 0.000 -0.333
domainCOG 0.000 -0.333 -0.333
dv_fctr_EXT -0.179 0.000 0.000 0.000
dmPOS:__EXT 0.000 -0.022 0.007 0.007 0.000
dmNEG:__EXT 0.000 0.007 -0.022 0.007 0.000 -0.333
dmCOG:__EXT 0.000 0.007 0.007 -0.022 0.000 -0.333 -0.333
(Note that domain was effect-coded and leaves out the bodily sensations domain here, but we could re-code it in all the different ways to get the contrast between bodily sensations and the grand mean as desired.)
I see three take-aways from these regression results, all of which I think I can see reflected in the plots above:
d_dev_factor_other %>%
mutate(domain = recode_factor(
domain,
"BOD" = "Bodily sensations",
"NEG" = "Negative emotions",
"POS" = "Social abilities & positive emotions",
"COG" = "Cognition & control",
.default = NA_character_),
capacity = gsub("_", " ", capacity)) %>%
filter(response != "idk") %>%
distinct(domain, capacity, response) %>%
arrange(domain, capacity, response) %>%
kable(caption = "Free responses to optional question about other developmental factors that might play a role (excluding answers like 'N/A,' 'no,' 'not that I can think of', etc.)") %>%
kable_styling() %>%
collapse_rows(1:2)
| domain | capacity | response |
|---|---|---|
| Bodily sensations | feeling pain | Again, you're born with this ability |
| I think feeling pain is innate | ||
| I think pain is one of the first sensations that humans are able to actually comprehend | ||
| It is just a natural sensation. Our body is programmed to feel bad to avoid further injury or notifiy us that one has occured | ||
| Nervous system | ||
| None come to mind | ||
| Not that I can think of. I think everyone is capable of feeling pain, even at birth, as it is something physical. Whether they understand why they feel it though is a different story | ||
| Pain is likely the very first thing we all experience when we are born. Childbirth is not exactly pleasant for the mother or the baby | ||
| This is another one I would use the word "primal" to describe | ||
| Watching other's reaction to painful situations | ||
| Whether or not they could ignore it or recall feeling pain before and cope | ||
| getting hungry | All children are capable of being hungry | |
| Biological Need | ||
| child is controlling himself | ||
| everyone gets hungry no matter what the age is . it depends on how we show we are hungry | ||
| evolution, survival mechanisms | ||
| Getting hungry is one of the first things a newborn will experience | ||
| I think we're all just born with the feeling of hunger because he have to have nourishment | ||
| It's just survival at this point | ||
| just something theyre born with | ||
| None come to mind | ||
| the child being active and using enagy | ||
| The child has emotions that cause them to increase their desire to do this | ||
| the child just gets hungry- human nature | ||
| the more experience they have with food, the bigger the capacity for hunger | ||
| We are just preprogrammed to know when we are hungry. | ||
| Yeah, I think this one could be covered with just one answer and needs no further thought or explanation | ||
| You're born with this ability. Technically, you have this ability BEFORE you're born | ||
| Negative emotions | feeling distressed | again this is similar to the first one as they age they become more independent and spend less time distressed |
| Babies show that they can be distressed from birth, or even during labor. | ||
| Distress is a biological survival mechanism | ||
| Distress is easily felt regardless of age. The factors change. A newborn gets distressed if they are hungry, or tired, or sick, etc. A 5 year old can get distressed if their favorite toy breaks | ||
| Feeling distressed is natural as the child has physical needs and cannot attend to them on their own | ||
| General living | ||
| I believe the feeling of distress is another innate ability that we are able to feel from the moment we're born. I think living, breathing creature that has feelings is the same way from birht | ||
| I do not imagine that any child of any age does not have this ability | ||
| I feel like feeling distressed can be attributed to human's evolution. Distress would be hard coded in the DNA (panic when something bad happens help survival) | ||
| I think babies crying shows distress | ||
| I think it's similar to the pain feeling. They can feel distressed, but they may not understand the feeling | ||
| I think that feeling distressed is somewhat primal. Flight or flight, so to speak. I don't know if I would use the word "biological" to describe it, given my choice, but it's the closest out of the option provided. Then there are layers of socialization, etc. beyond that which amplify it | ||
| i think they always know distress in my opinion | ||
| if the child is left alone or isnt care for | ||
| It's an automatic response | ||
| It's the result of biological function, like pain response | ||
| Maybe hereditary? | ||
| Maybe the amount of nutrition a child gets throughout their life cycle | ||
| None come to mind | ||
| There could be physiological elements that prevent a child from being able to feel distress as acutely as other children | ||
| They are alive | ||
| They can feel stress from the people around them especially parents even inside the womb | ||
| You're born with all emotions, even if you don't understand them yet | ||
| feeling helpless | Again, you're born with emotions. You might not fully comprehend them, but you have the ability to feel them, nontheless | |
| as they age there is a tendency to become more independent and at that age less aware of danger | ||
| Feeling helpless is something that is learned through experience, not an innate ability | ||
| Feeling helplessness is not really an ability at all, it's an emotion. Most people are capable of feeling helpless | ||
| I don't think you can feel helpless until you become completely self aware, and I'm not sure when that happens | ||
| I feel like a kid would need to have some understanding of their emotions | ||
| I feel that it changes. A younger child may not know or understand the feeling of helpless. An older child my feel they are invinsible or independent. | ||
| I think helplessness stems from a fundamental fear of losing something we've got or not getting something we want (regardless of how well we are able to articulate or understand what that is in a sophisticated way) | ||
| I think it is just mental development. Young children don't really feel self pity | ||
| I think that it all surrounds object perception/realization. Once the child realizes that "food" makes them not hungry or "moving" makes them not hurt | ||
| I think this is something that happens later in life | ||
| My zero is because at the beginning, there's no feeling of self, so you can't feel helpless if you don't feel the self (in my opinion) | ||
| None come to mind | ||
| Not being taken care of or given what it wants would eventually make them real familiar with feeling helpless | ||
| The child has negative experiences that increase this | ||
| The child notices when it is not cared for properly by feelings of hunger,pain etc | ||
| The child's need to be loved and cared for contributes to their feeling helpless, since they are generally moreso helpless in the earlier stages of life than later | ||
| The child's needs change | ||
| This is also hard to rate but I go with my same rationale for happiness. Infants ARE helpless and I think on some leel they know it | ||
| Social abilities & positive emotions | feeling happy | A child that feels love feels happiness |
| Again, I think happiness (albeit, not amusement) is something relatively innate related to security & safety. I feel like this should be an incredibly easy question to pinpoint to a specific cause, and yet somehow it is not at all | ||
| As the brain grows, we are capable of sensating new emotions. | ||
| I think it's innate in children..Sometimes they interpret it differently as they get older though | ||
| I think its just within them to feel happy | ||
| it is a normal emotion | ||
| None come to mind | ||
| Other external stimuli encourages this | ||
| social environment | ||
| The child becomes more able to notice their own feelings | ||
| the child senses the love shared around him\her | ||
| The child's interests/needs change | ||
| This is a harder one for me to rate. I would say some form of emotional range is present at birth. It would have to be. Beyond that I do not know | ||
| Well fed or not | ||
| yes it depends if the child is in a good home or a bad home | ||
| You're born with this ability | ||
| learning from other people | after the age of 2 or 3 the ability to learn langauges begins to diminish hence my ratings after that age | |
| I feel like this is just the nature of human beings to learn from the others around (we wouldn't be where we are without it) | ||
| I think there are a ton of factors that lend to learning. I think it was one of the major purposes of the brain and so many, many parts of the brain and body contribute to it | ||
| kids love to copy people they are around and that how they learn | ||
| Learning is an innate ability. Generally, everyone is born with the ability to learn | ||
| None come to mind | ||
| Relationship with parents | ||
| The amount of nutrition the child gets throughout their life cycle | ||
| The more they age, the more connections the brain makes. | ||
| Cognition & control | controlling their emotions | A 0 to 9 month old infant is never wrong. You don't tell a 0 to 9 month old infant "No" or "Not right now". After 9 months, then you can start having that conversation |
| Children imitate other children | ||
| having very good role models and a good home to learn these in | ||
| I don't think anyone can control thier emotions (make yourself sad, NOW!), only the expression of them | ||
| I know a lot of adults that can't control their emotions. This is something you either have or don't have, | ||
| I would say the type of parenting being done by the parents | ||
| I'm afraid I don't have much to add here. Most of your questions covered it. I would say their environment, but that was, more or less, covered in the questions | ||
| Learning from adults reaction to behavior | ||
| Learning that outcomes can change depending on the the emotion they display | ||
| None come to mind | ||
| overall social environment | ||
| The child has life experiences after birth that increase their ability to do this | ||
| They aren't necessarily learning to control the emotion so much as the expression | ||
| They brains haven't fully developed they cab;'t help themselves and are still exploring their emotions | ||
| They react to positive or negative responses from other humans, if they control, or don't control their emotions | ||
| watching tv or movies | ||
| reasoning about things | Again, up to about 9 months old the wiring is just not there. It beginbs to develop but... to be honest reasoning is something many adults fail at | |
| Genetic driven intelligence ability | ||
| injury/disability | ||
| It's mostly just natural mental development | ||
| Maybe the amount of food/nutrition they get throughout their life cycle | ||
| None come to mind | ||
| not sure if they can have those skills | ||
| Nuerological development can impact this..for example if they develop brain injuries...autism etc | ||
| The more they age, the more they comprehend and understand. | ||
| watching tv or videos |
ggplot(d_dev_factor_most_important_choice %>%
mutate(domain = recode_factor(
domain,
"BOD" = "Bodily sensations",
"NEG" = "Negative emotions",
"POS" = "Social abilities & positive emotions",
"COG" = "Cognition & control",
.default = NA_character_)) %>%
mutate(dev_factor = case_when(
grepl("Other", response) ~ "other",
grepl("teach", response) ~ "people teach",
grepl("experiments", response) ~ "experiments",
grepl("womb", response) ~ "womb experiences",
grepl("interacts", response) ~ "interacts people",
grepl("preprogrammed", response) ~ "preprogrammed",
grepl("objects", response) ~ "observes objects",
grepl("observes the people", response) ~ "observes people",
grepl("body grows", response) ~ "body grows",
grepl("brain changes", response) ~ "brain changes",
grepl("senses improve", response) ~ "senses improve"),
dev_factor = factor(
dev_factor,
levels = c(gsub("_", " ",
levels(d_dev_factor_rating$dev_factor)),
"other"))) %>%
mutate_at(vars(capacity, response),
funs(gsub("_", " ", .))),
aes(x = dev_factor, fill = domain)) +
facet_wrap(~ domain ~ capacity, ncol = 4) +
theme(legend.position = "none",
axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
geom_bar() +
labs(title = "Choice of most important developmental factor (Version 1)",
x = "Developmental factor",
y = "Number of participants")
temp_other <- d_dev_factor_most_important_choice %>%
filter(grepl("Other", response)) %>%
count(domain, capacity)
ggplot(d_dev_factor_most_important_choice %>%
mutate(domain = recode_factor(
domain,
"BOD" = "Bodily sensations",
"NEG" = "Negative emotions",
"POS" = "Social abilities & positive emotions",
"COG" = "Cognition & control",
.default = NA_character_)) %>%
mutate(dev_factor = case_when(
grepl("Other", response) ~ "other",
grepl("teach", response) ~ "people teach",
grepl("experiments", response) ~ "experiments",
grepl("womb", response) ~ "womb experiences",
grepl("interacts", response) ~ "interacts people",
grepl("preprogrammed", response) ~ "preprogrammed",
grepl("objects", response) ~ "observes objects",
grepl("observes the people", response) ~ "observes people",
grepl("body grows", response) ~ "body grows",
grepl("brain changes", response) ~ "brain changes",
grepl("senses improve", response) ~ "senses improve"),
dev_factor = factor(
dev_factor,
levels = c(gsub("_", " ",
levels(d_dev_factor_rating$dev_factor)),
"other"))) %>%
mutate_at(vars(capacity, response),
funs(gsub("_", " ", .))) %>%
filter(dev_factor != "other"),
aes(x = reorder(capacity, as.numeric(domain)), fill = domain)) +
facet_wrap(~ dev_factor, ncol = 5) +
theme(legend.position = "top",
axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
geom_bar() +
labs(title = "Choice of most important developmental factor (Version 2)",
subtitle = paste0("Excluding answers of 'other' (n < ",
max(temp_other$n) + 1, ") for any capacity"),
color = "Domain",
x = "Capacity (by domain)",
y = "Number of participants")
d_dev_factor_most_important_free %>%
mutate(domain = recode_factor(
domain,
"BOD" = "Bodily sensations",
"NEG" = "Negative emotions",
"POS" = "Social abilities & positive emotions",
"COG" = "Cognition & control",
.default = NA_character_),
capacity = gsub("_", " ", capacity)) %>%
distinct(domain, capacity, response) %>%
arrange(domain, capacity, response) %>%
kable(caption = "Free responses when 'other' was selected as most important developmental factor") %>%
kable_styling() %>%
collapse_rows(1:2)
| domain | capacity | response |
|---|---|---|
| Bodily sensations | feeling pain | everyone is capable of feeling pain. |
| getting hungry | All children are capable of being hungry. | |
| Negative emotions | feeling distressed | As soon as his consious develops they are considered alive and capable of feeling. |
| emotional development | ||
| most id not all children are capable of experinceing distress not matter what. | ||
| feeling helpless | Experiences in life | |
| I really don't know. | ||
| Not sure helplessness is felt so young | ||
| Social abilities & positive emotions | feeling happy | Every child is capable of feeling happiness no matter what. |
| it is just a normal response too something enjoyed | ||
| learning from other people | not sure | |
| Cognition & control | controlling their emotions | impossibe |
ggplot(d_demo, aes(x = Duration/60)) +
geom_histogram(binwidth = 2) +
geom_vline(xintercept = median(d_demo$Duration/60), color = "blue", lty = 2) +
scale_x_continuous(breaks = seq(0, 10000, 4)) +
labs(title = "Duration of study (according to Qualtrics)",
subtitle = "Blue dotted line marks median",
x = "Duration (in minutes)",
y = "Number of participants")
ggplot(d_demo, aes(x = Age)) +
geom_histogram(binwidth = 2) +
geom_vline(xintercept = median(d_demo$Age), color = "blue", lty = 2) +
scale_x_continuous(breaks = seq(0, 10000, 4)) +
labs(title = "Particpiant age (self-reported)",
subtitle = "Blue dotted line marks median",
x = "Age (in years)",
y = "Number of participants")
ggplot(d_demo, aes(x = GenderSex)) +
geom_bar() +
labs(title = "Particpiant gender/sex (self-reported)",
x = "Gender/sex",
y = "Number of participants")
ggplot(d_demo, aes(x = gsub('(.{1,30})(\\s|$)', '\\1\n', RaceEthnicity_collapse))) +
geom_bar() +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
labs(title = "Particpiant race/ethnicity (self-reported)",
x = "Race/ethnicity",
y = "Number of participants")
ggplot(d_demo, aes(x = FirstLang)) +
geom_bar() +
labs(title = "Particpiant first language (self-reported)",
x = "Language",
y = "Number of participants")
ggplot(d_demo, aes(x = factor(Education,
levels = levels(d$Education),
labels = gsub('(.{1,30})(\\s|$)', '\\1\n',
levels(d$Education))))) +
geom_bar() +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
labs(title = "Particpiant educational attainment (self-reported)",
x = "Highest level of education completed",
y = "Number of participants")
ggplot(d_demo, aes(x = Income)) +
geom_bar() +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
labs(title = "Particpiant household income (self-reported)",
x = "Annual household income",
y = "Number of participants")
ggplot(d_demo, aes(x = HouseholdSize)) +
geom_histogram(binwidth = 1) +
geom_vline(xintercept = median(d_demo$HouseholdSize), color = "blue", lty = 2) +
scale_x_continuous(breaks = seq(0, 10000, 1)) +
labs(title = "Particpiant household size (self-reported)",
subtitle = "Blue dotted line marks median",
x = "Number of people in household (adults and children)",
y = "Number of participants")
ggplot(d_demo, aes(x = MaritalStatus)) +
geom_bar() +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
labs(title = "Particpiant marital status (self-reported)",
x = "Marital status",
y = "Number of participants")
ggplot(d_demo, aes(x = Parent)) +
geom_bar() +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
labs(title = "Particpiant parent status (self-reported)",
subtitle = "'NA' indicates response of 'Prefer not to say'",
x = "Parent status",
y = "Number of participants")
ggplot(d_demo %>% filter(Parent == "Yes"), aes(x = ChildrenNumber)) +
geom_histogram(binwidth = 1) +
geom_vline(xintercept = median(d_demo[d_demo$Parent == "Yes",]$ChildrenNumber, na.rm = T),
color = "blue", lty = 2) +
scale_x_continuous(breaks = seq(0, 10000, 1)) +
labs(title = "Number of children among parents (self-reported)",
subtitle = "Blue dotted line marks median",
x = "Number of children (among parents)",
y = "Number of participants")
ggplot(d_demo %>% filter(Parent == "Yes"),
aes(x = factor(ChildrenYoungestAge_collapse,
levels = levels(d_demo$ChildrenYoungestAge_collapse),
labels = gsub('(.{1,30})(\\s|$)', '\\1\n',
levels(d_demo$ChildrenYoungestAge_collapse))))) +
geom_bar() +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
labs(title = "Age of youngest child among parents (self-reported)",
x = "Age of child",
y = "Number of participants")
ggplot(d_demo %>% filter(Parent == "Yes"),
aes(x = factor(ChildrenOldestAge_collapse,
levels = levels(d_demo$ChildrenOldestAge_collapse),
labels = gsub('(.{1,30})(\\s|$)', '\\1\n',
levels(d_demo$ChildrenOldestAge_collapse))))) +
geom_bar() +
theme(axis.text.x = element_text(angle = 45, hjust = 1, vjust = 1)) +
labs(title = "Age of oldest child among parents (self-reported)",
x = "Age of child",
y = "Number of participants")
Social abilities & positive emotions
Linear effects only
Non-linear effects